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Se  /cra]  examples  of  the  numerical  evaluation  of  an  integral 
equation  for  the  calculation  of  the  attenuation  of  a  radio  wave  are  given. 
These  waves  are  assumed  to  be  propagated  over  realistic,  smoothly 
varying  irregular,  inhomogeneous  terrain.  Results  for  propagation 
over  a  cylindrical  earth  show  an  accuracy  to  3-4  significant  figures 
when  compared  with  the  classical  residue  series.  Calculations  for 
propagation  over  smooth  mixed  leUid-sea  paths  agree  with  classical 
methods.  The  applicability  of  the  program  to  permit  computation  of 
propagation  over  terrain  with  smooth  height  variation  is  demonstrated 
by  calculations  of  propagation  over  one  and  two  GaussicUi-shaped  hills. 
The  ability  of  the  program  to  allow  treatment  of  variations  in  both 
ground  cotiductivity  and  height  combirod  is  illustrated  by  calculations 
of  propagation  from  the  sea  up  a  sloping  beach  and  by  calculations  of 
propagation  over  an  island.  This  last  example  illustrates  the 
importance  of  the  terrain  profile  in  mixed  path  calculations. 
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1.  INTRODUCTION 


Despite  numerous  attempts,  a  numerically  feasible  way  to 
calculate  the  field  strength  of  a  radio  wave  propagating  over  realistic, 
smoothly  varying,  inhomogeneous  terrain,  has  not  been  found. 
Hufford  (1952)  developed  an  integral  equation  for  such  propagation  by 
using  the  free-space  Green's  function  in  Green's  second  identity  and 
showed  that  his  solution  yielded  the  c’assical  result  for  propagation 
over  a  smooth  sphere.  Berry  (1967)  succeeded  in  solving  the  equation 
numerically  for  vertically  polarized  radio  waves,  showing  sample  cal¬ 
culations  up  to  10  MHz.  If  the  normalized  surface  impedance  is  not 
much  smaller  than  1,  the  numerical  techniques  are  very  inefficient, 
however,  and  round-off  errors  accumulate  so  fast  that  the  results  are 
not  useful.  For  normal  ground  constants,  this  condition  excludes  all 
horizontally  polarized  waves  and  all  vertically  polar irrd  waves  above 
a  few  megahertz. 

The  method  used  in  this  paper  is  based  on  an  elementary  function 
that  is  closely  related  to  the  Sommerfeld  flat  earth  attenuation  function. 
This  elementary  function  satisfies  a  scalar  'parabolic''  wave  equation. 
The  resulting  integral  equation  is  numerically  feasible  for  both  vertical 
and  horizontal  nolarization  and  for  normalized  surface  impedances  in 
the  HF  band. 

The  problem  to  be  solved  is  illustrated  in  figure  1,  which  shows  a 
possible  propagation  path.  The  signal  at  the  receiver  is  affected  by  the 
mean  curvature  of  the  earth,  height  profile  along  the  path,  and  the  ch  -ng 
of  surface  impedance  along  the  path.  The  changes  may  be  abrupt  (o.g., 
at  a  land  sea  boundary),  or  gradual  (e.g.,  as  the  sea  state,  temperature 
or  salinity  change).  The  problem  of  abrupt  changes  in  surface  imped¬ 
ance  at  smooth  land-sea  boundaries  has  been  solved  (Wait, 

Numerical  results  for  changes  in  surface  impedance  have  been  calcu¬ 
lated  by  Rosich  (1968,  1970). 


The  present  work  allows  the  terrain  to  be  represented  by  a 
completely  arbitrary  profile  in  terms  of  the  elevation  versus  distance. 
The  hills  and  valleys  themselves  are  taken  to  be  uniform  in  the  direction 
transverse  to  the  propagation  direction.  The  terrain  may  also  be  char¬ 
acterized  by  a  conductivity  and  dielectric  constant  which  are  functions  of 
distance. 

The  main  body  of  the  report  describes  the  results  of  calculations 
for  several  examples  including  paths  similar  to  that  in  figure  : .  The 
appendices  contain  the  derivation  of  the  integral  equation,  the  necessary 
numerical  analysis,  and  a  listing  of  the  Fortran  computer  program. 


2.  THE  INTEGRAL  EQUATION 


The  derivation  of  the  integral  equation  is  given  in  appendix  A.  The 
details  will  not  be  reiterated  here,  but  the  final  result  is  (Otl,  1971  ) 

X 


f(x)  =  g(x,y)  W(x,  0)  Iy  J  jyM?)  W(x,  ?)  - 

V 


y(x)  -  y(^ 

X  • 


+  (a  -  4.) 


r  /  L  '  (X-  ?)  J 


d«r  , 


(1) 


where  x,  »,  y(x)  and  y(»)  are  defined  in  figure  The  factor  (  •  -  \  ^) 
arises  in  mixed-path  problems.  That  is,  substituting  for  i  in  |A-2» 
and  (A-9)  will  yield  the  difference  (^-  The  factor  ^  is  c«.inslant 

with  distance  and  is  the  relative  value  of  the  normalized  surface 
impedance.  This  factor  is  computed  using  the  values  for  -  and  ^  for 
the  first  section  of  a  mixed  path.  The  factor  •  varies  with  distance  in 
a  mixed  path  problem.  The  variation  of  >  with  x  may  be  continuous  <>r 
contain  abrupt  changes.  The  factor  (^-  *,  )  is  zero  tor  a  single  section 

path.  The  remaining  factors  in  (1)  are  defined  as 
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U)(X.  ?)  = 


$ 


[y(3c)  -  y(g)f  +  ifili  _  y^W 

2(x  -  ^  2e  2x 

W(X.  5  )  =  1  -  i  /np  w  (-  /TT  )  , 


p  =  -ik  /^»(x  -  5)/2  , 
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w(-  ../TT )  =  e  erfc  (i  /VT ) 


1  •  "**  dt 

=  ~  J  ^  (Abramowitr  and  Stegun,  1964) 


&  = 


vertical  polarisation 
horizontal  polarization 


il8(10»)o 
^  ®r  £(MHz) 

f  s  frequency,  in  MHz 

(7  -  ground  conductivity 

e  -  dielectric  constant 

g(^  y)  =  antenna  pattern  factor. 


Equation  (1)  gives  the  integral  equation  for  the  attenuation  function 
normalized  to  tw  ice  the  free-space  field.  The  details  of  the 
numerical  solution  of  (1)  are  given  in  appendix  B.  Since  the  upper 
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limit  of  integration  in  (1)  is  x,  the  effects  of  backscatter  are  excluded. 
That  is,  to  include  the  effects  of  backscatter,  the  range  of  integration 
would  have  to  include  the  entire  terrain.  Also,  the  integral  equation 
in  (1)  neglects  the  effects  of  side-scatter  since  the  derivation  of  (1) 
assumed  ridges  uniform  in  the  direction  transverse  to  the  propagation 
direction.  In  the  case  of  small  slopes  and  the  transmitting  antenna 
near  the  earth,  side-scatter  and  backscatter  are  second  order  effects. 

3.  EXAMPLES 

In  this  section  we  examine  the  behavior  of  the  attenuation  function, 
f(x), for  eight  terrain  profiles,  y(x).  Comparisons  of  results  from  (1) 
with  previous  results  for  a  flat  earth,  a  smooth  homogeneous  cylindrical 
earth,  a  smooth  sea-land- sea  path  and  a  single  Gaussian- shaped  ridge 
seem  to  validate  the  technique.  Its  more  general  applicability  is  illus¬ 
trated  by  calculations  for  propagation  over  two  Gaussian  hills,  over  an 
island  that  rises  above  sea  level,  and  over  a  sea- sloping  beach  with  a 
Sana- dune  path, 

3. 1  A  Flat  Earth 

y(x)  =  0,  y'(x)  =  0,  The  solution  of  the  integral  equation  (1)  is 
trivial  and  is  simply 

f(x)  =  W(x),  (2) 

where  W(x)  is  the  Sommerfeld  flat- earth  attenuation  function  (Wait, 
1964). 

3.2  A  Parabodoidal  Earth 

y(x)  2-  -x3 /2a,  y'(x)  =  -x/a,  where  a  is  the  radius  of  the 
cylinder  and  is  taken  to  be  about  6.  37  x  10^  kilometers.  The  frequency 
of  the  transmitting  antenna  is  1  MHz  and  is  vertically  polarized.  The 
ground  constants  are:  (j  =  0.01  mho/m  and  0^  =  10.  The  magnitude 
and  phase  of  the  attenuation  function  versus  horizontal  distance,  x  are 
given  in  table  I,  These  results  are  compared  in  table  I  with  those 
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obtained  using  the  residue  series  (V/ait,  1964)  for  the  attenuation  function. 
The  agreement  is  seen  to  be  very  good  out  to  the  largest  distance  com¬ 
puted.  F'  r  example,  at  300  km,  the  difference  in  the  phase  between  the 
integral  equation  method  and  the  residue  series  method  is  about  0.  009 
rad.  or  about  0.  5“.  The  greatest  error  in  amplitude  occurs  at  a  distance 
of  150  km  and  is  about  2  units  in  the  third  significant  figure.  The  error 
decreases  on  either  side  of  this  point,  a  characteristic  common  to  many 
numerical  solutions.  The  results  obtained  in  table  I  are  for  a  step  size, 
h  =  1  km;  hovrsver,  a  step  size  of  2  km  did  not  change  the  results 
appreciably.  A  detailed  error  analysis  is  beyond  the  scope  of  this  paper. 
The  last  significcint  figure  of  agreement  in  table  I  is  underlined. 

Table  I.  Atteniiation  function  versus  distance 


Integral  Equation  Soluticm  Residue  Series  Solution 


Horizontal 
Distance,  x 
(km. ) 

Amplitude 

Phase 
(rad. ) 

Amplitude 

Phase 
(rad. ) 

0 

1.0 

0 

1.0 

0 

25 

0. 51331 

-1.9717 

0. 51332 

-1.9709 

50 

0.28936 

-2. 5929 

0.28970 

-2. 5921 

75 

0.17575 

-2.9^97 

0. 17595 

-2.9556 

100 

0. 11506 

3.0902 

0.11520 

3.0892 

125 

0.  08£44 

2.9126 

0. 08000 

2.9131 

150 

0.  05914 

2.7^06 

0.05939 

2.  7663 

175 

0. 04504 

2.  6J.69 

0.04502 

2.  6120 

200 

0. 035^2 

2.4736 

0.03509 

2.4680 

225 

0.  027811 

2.  3278 

0.02777 

2.  3213 

250 

0.022^4 

2. 1J80 

0.  02221 

2. 1710 

275 

0.01790 

2.  0249 

0.  01788 

2.  01 68 

300 

0.01442 

1.8281 

0.01446 

1.8591 
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The  time  required  to  compute  the  attenuation  function  at  inter¬ 
vals  of  one  kilometer  out  to  a  maximum  distance  of  300  kin  was  about 
25  min  using  a  CDC  3800  computer.  The  time  required  to  compute  the 
attenuation  function  for  a  specified  profile  is  approximately  proportional 
to  the  square  of  the  number  of  points  used  along  the  abscissa.  Thus,  in 
the  above  example,  if  the  maximum  distance  were  1 50  km  rather  than 
300  km,  the  time  required  would  be  about  1/4  as  much,  or  about  6  min. 
The  sample  input  and  output  data  given  in  appendix  C  pertain  to  this 
sample. 

3.  3  A  Gaussian -Shaped  Ridge 

-(x-5f 

y  =  e  t  y*  -  -2{x-5)y.  This  is  a  more  interesting  profile 

at  least  from  the  standpoint  of  radio  propagation.  The  profile  together 
with  the  magnitude  of  the  attenuation  function  versus  distance  are  shown 
plotted  in  figure  3.  The  magnitude  of  the  attenuation  function  |  f(x)|  , 
is  normalized  to  twice  the  free  space  field,  2  exp( -ikr^)  /  r^.  The  ob¬ 
server  is  located  on  the  terrain  and  the  transmitter  is  located  at  the 
coordinate  origin.  The  ground  constants  are  a  =  0.  01  mho/  m  and 

e  a  10,  The  transmitter  is  vertically  polarized  zmd  Ihe  frequency  is 
r 

1  MHz.  The  terrain  profile  shown  in  the  insert  has  a  maximum  height 
of  1  km  and  the  hill  is  centered  at  a  point  5  km  from  the  transmitter. 
The  solid  straight  line  in  figure  3  is  the  attenuation  function,  W(x), 
for  a  flat  earth. 

The  data  in  figure  3  represented  by  crosses  was  obtained  by 
replacing  the  Gaussian-shaped  ridge  with  a  rounded  knife-edge 
and  computing  the  field  on  the  surface  shown  dashed  in  figure  4 
using  "4-ray  theory"  (Schelleng,  et.  al. ,  1933).  The  radius  of  the 
rounded  knife-edge  is  500  m  (which  is  the  curvature  of  the  Gaussizm 
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Hufford 

itegral  Equation  Sommerfeld  Flat  Earth 


(receiving  antenna)  is  located  on  the  profile. 
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hill  at  its  crest)  and  the  knife-edge  is  located  1  km  above  the  plane 
y  =  0.  The  four  rays  are  the  two  rays  that  strike  the  knife-edge  on  the 
illuminated  side  plus  the  two  rays  that  reach  an.  observer  in  the  shaded 
side  i.  e. ,  a  direct  diffracted  ray  and  a  ray  which  is  diffracted  euid  then 
reflected  before  reaching  the  observer.  The  results  in  figure  3  show 
excellent  agreement  between  the  points  computed  using  "4-ray  theory" 
and  those  obtained  solving  the  integral  eqtiation  numerically. 

The  open  circles  in  figure  3  were  computed  using  the  Hufford 
integral  equation  (Hufford,  1952).  Since  there  aro  fewer  approximations 
in  the  Hufford  integral  equation  than  in  the  results  presented  in  thic 
paper,  the  former  shovild  be  considered  the  most  accurate.  Hufford*  s 
integral  equation  shows  a  slight  dip  in  the  attenuation  function  at  a 
distance  of  about  9  km  which  is  exaggerated  by  the  solid  circles  but 
does  not  appear  in  the  knife-edge  results.  Also,  the  open  circles  dif¬ 
fer  somewhat  in  the  shadow  from  the  results  presented  earlier  by 
Berry  (1967).  There  were  projection  factors,  of  the  form  >J  I  +  {y*f  , 
omitted  from  Berry’ s  results  since  in  most  applications  these  factors 
are  nearly  unity,  i.e,,  y'  is  small.  However,  in  the  present  example 
these  factors  become  important. 

The  solid  circles  in  figure  3  present  the  attenuation  function 
computed  numerically  using  the  integral  equation  in  (1).  We  find  some 
error  in  the  results  obtained  using  the  integral  equation  presented  in 
this  paper  around  6  km  and  9  km.  However,  the  error  is  small  and  is 
exciggerated  in  this  particiilar  example  because  of  the  large  slopes 
encountered  on  the  terrain  profile.  The  error  is  a  result  of  the 
as  sumption  that 

a  0  , 

or  that  the  fast  phase  variation  of  cp  with  x  is  in  the  term  exp(-ikx). 
In  most  terrain  profiles,  this  will  indeed  be  a  good  approximation  and, 
in  fact,  in  the  present  example  yields  adequate  accuracy. 


d  X® 
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The  physical  characteristics  of  the  results  in  figure  3  are  inter¬ 
preted  most  easily  using  the  ray  picture.  The  attenuation  function 
decreases  at  the  flat  earth  rate  for  the  first  km.  Then*  as  the 
observation  point  moves  up  the  crest  of  the  hill,  the  attenuation  function 
increases  due  to  focusing  of  the  direct  ray  and  the  surface  ray  on  the 
lit  side  of  the  crest.  The  attenuation  function  reaches  its  maximum 
value  very  close  to  the  point  on  the  terrain  uhere  there  is  an  inflection 
point.  This  increase  in  the  amplitude  to  a  maximum  on  the  lit  side  near 
the  crest  has  also  been  predicted  analytically  by  Wait  and  Murphy 
(1938).  Just  over  the  top  of  the  hill  the  attenuation  function  decreases 
since  the  direct  ray  is  no  longer  present  and  then  the  attenuation  func¬ 
tion  partially  recovers  again  due  to  the  constructive  interference  of  a 
direct  diffracted  ray  and  a  diffracted  ray  traveling  along  the  surface 
before  reaching  the  observer. 

3. 4  A  Sea  -Land-Sea  Path 

The  terrain  profile  is  fiat  in  this  example  and  the  ground  con¬ 
stants  change  abruptly  at  the  sea-land,  land  sea  interfaces.  This 
example  was  selected  as  a  check  on  the  mixed  path  capabilities  of  the 
method.  The  results  for  the  magnitude  of  the  attenuation  function 
normalized  to  twice  the  free  space  field  are  plotted  in  figure  5  versus 
distance  from  the  antenna  in  km.  The  antenna  is  vertically  polarized 
and  the  frequency  is  10  MHz.  The  solid  circles  represent  the  attenua¬ 
tion  function  computed  numerically  using  (1).  The  open  circles  in 
figure  5  represent  the  attenuation  function  computed  by  Rosich  (1968, 
1970)  using  a  perturbation  approach.  The  data  given  by  the  crosses  in 
figure  5  represents  the  attenuation  function  computed  using  a  method 
based  upon  &e  classical  residue  series  (Furutsu,  et.  al. ,  1964).  This 
method  is  eq\iivalent  to  that  of  Wait  (1964).  This  latter  method  makes 
the  fewest  approximations  for  the  three  section  earth  considered  in 
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|f  (x)|,  Mognitude  of  Attenuation  Function  Normalized  to  2 exp  (-ikro)/r 


X,  km 

Figur«  5.  A  path.  Tha  profila  li  flat. 
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this  example.  The  agreement  between  the  solid  circles  representing  (1) 
and  the  crosses,  appears  to  demonstrate  the  validity  of  the  formulation 
in  treating  mixed  path  propa^tion  problems.  The  abrupt  changes  in 
conductivity  and  dielectric  constant  used  in  this  example  do  not  repre¬ 
sent  a  realistic  sea-land  interface.  The  method  is,  however,  capable  of 
treating  a  continuous  variation  of  conductivity  and  dielectric  constant. 

3.  5  A  Sea-Land-Sea  Path  With  An  Island 

This  example  combines  terrain  features  and  mixed-path  effects. 
The  island  is  drawn  to  scale  in  figure  6  and  its  elevation  is  250  m  at  the 
highest  point.  The  magnitude  of  the  attenuation  function  normalized  to 
twice  the  free  space  field  versus  distance  is  plotted  in  figure  6.  The 
antenna  is  vertically  polarized  and  the  frequency  is  10  MHz.  For  com¬ 
parison,  the  magnitude  of  ttic  attenuation  function  for  a  flat  island  is 
also  shown  in  figure  6.  The  most  significant  feature  of  figure  6  is  that 
the  terrain  profile  has  a  greater  effect  on  the  attenuation  function  on  the 
island  than  do  changes  in  the  ground  constants,  and  the  residual  effect 
of  the  profile  well  beyond  the  island  is  comparable  to  that  of  ttxe  change 
in  ground  constants. 

3.  6  A  Sloping  Beach  At  High  And  Low  Tides 

The  profile  is  drawn  to  scale  in  figure  7  and  the  assumed  ground 
constants  u)ed  for  the  wet  and  dry  sand  are  given  in  Qie  figure.  The 
transmitter  is  out  at  sea.  As  the  tide  rises,  the  wet  sand  in  figure  7  is 
covered  by  water  and  as  the  tide  recedes  it  exposes  the  wet  sand.  The 
magnitude  of  the  attenuation  function  versus  distance  is  shown  plotted  in 
figure  7.  There  is  little  difference  in  the  attenuation  functioi  at  high 
and  low  tide.  However,  the  presence  of  the  crest  in  ths  beach  produces 
a  peak  in  the  atteiiuation  function  on  Uie  lit  side  and  a  shadow  in  back. 
This  illustrates  the  importance  of  the  terrain  profile  in  mixed  path 
problems. 
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Figure  7.  A  eloping  bench  at  high  and  low  tidee, 


3. 7  Two  Gaussian  Hills 


The  profile  is  drawn  to  scale  in  figure  8*  The  separation  of  the 
hills  is  such  that  a  null  instead  of  a  peak  in  the  attenuation  function  is 
produced  on  the  lit  side  of  the  second  hill.  Obviously  there  are  an  infinite 
number  of  combinations  of  hills  that  will  in  turn  produce  an  infinite 
number  of  possible  combinations  of  nulls  and  peaks  in  the  attenviation 
function.  The  method  will,  in  principle,  treat  any  nximber  of  hills  and 
valleys.  The  hills  need  not  have  Gaussian  profiles;  any  smooth  func¬ 
tion  of  distance  is  acceptable. 

3. 8  A  Gaussian  Hill  (transmitting  frequency  of  10  MHz  ) 

The  profile  as  well  as  the  magnitude  of  the  attenuation  fimction 
versus  distance  is  shown  in  figure  9.  The  results  in  figure  9  differ 
somewhat  from  those  published  earlier  by  Berry  (1967).  Near  the 
crest  of  the  hill  small  oscillations  in  the  attenuation  occur  which  were 
not  present  when  the  transmitting  frequency  was  1  MHz.  One  possible 
explanation  for  these  wiggles  is  numerical  instability.  However,  this 
explanation  was  discarded  when  finer  subdivisions  of  the  integration 
interval  failed  to  remove  the  oscillations.  At  present,  they  can  only  be 
explained  in  terms  of  an  interference  efiect  between  a  ground-reflected 
wave  and  the  ground  wave  the  former  being  stronger  at  10  MHz  than 
at  1  MHa.  This  case  represents  a  quasi  upper  limit  in  the  capability  of 
the  computer  program  in  terms  of  frequency  and  slopes.  That  is, 
higher  frequencies  can  be  treated  but  the  terrain  cannot  change  as  fast 
as  It  does  in  figure  9.  Conversely,  more  rapid  changes  in  terrain  can 
be  treated  provided  the  frequency  is  less  than  it  is  in  figure  9.  Suice 
the  slopes  in  figure  9  are  near  unity,  we  have  a  heuristic  uncertainty 
principle  for  our  computer  program 

y'  f  i  10  (MHz)  . 
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4.  RECOMMENDATIONS  AND  CONCLUSIONS 


The  numerical  evaluation  of  an  integral  equation  for  ttie  propaga¬ 
tion  of  radio  waves  over  irregular*  inhomogeneous  terrain  is  demon¬ 
strated  for  several  examples.  Some  of  the  examples  provide  a  realistic 
picture  of  the  attenviation  of  a  radio  'Arave  vdien  it  encounters  a  terrain 
anomaly,  such  as  a  large  conducting  ridge.  The  Gaussian -Hill  example 
at  1  MHz  yields  physical  insight  into  a  focusing  phenomenon  of  the  field 
just  before  the  crest  of  a  hill  &at  cannot  be  predicted  on  the  basis  of 
simple  diffraction  theory,  but  is  in  fact  predicted  by  the  numerical 
solution  of  the  integral  equation.  However,  ray  theory  in  a  concave 
region  with  multiple  reflections  may  work. 

It  appears  that  the  results  discussed  in  this  report  represent  a 
useful  tool  for  analyzing  the  attenuation  loss  of  a  radio  wave  as  it  en¬ 
counters  terrain  anomalies  such  as  hills,  valleys,  land-sea  interfaces, 
etc.  The  computer  program  for  this  analysis  is  listed  in  appendix  C. 
However,  there  are  improvements  that  should  be  stu*iied.  They  are 
listed  below  in  an  order  not  necessarily  representing  their  relative 
importance. 

1)  A  three-dimensional  model  of  the  terrain.  It  should  be 
determined  if  the  energy  follows  a  geodesic  and  if  the 
effects  of  transverse  curvature  are  important  or  not. 

2)  Since  the  solution  represented  by  the  integral  equations 
does  in  fact  represent  a  solution  of  the  wave  equation 
plus  boundary  conditions,  it  applies  to  VHF  frequencies 
as  well  as  HF  frequencies.  Consequently,  numerical 
techniques  should  be  studied  so  that  the  program  will 
handle  VHF  frequencies  efficiently. 
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3) 


Real  antennas  rather  than  an  idealized  point  source  with 
an  arbitrary  pattern  factor  should  be  investigated; 
especially  when  a  large  diflracting  obstacle  is  within 
the  first  Fresnel  zone  of  the  antenna. 
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APPENDIX  A:  Derivation  of  the  integral  equation 
Consider  a  solution,  cp,  of  the  wave  equation 

^  k^cp  = -2TTT(x,y)  ,  y>y(x) 

which  satisfies  an  impedance  boundary  condition  of  the  form 


(ii) 


5cp  _  ikA  cp 


y  =  y(x) 


where  cp  represents  the  vertical  component  of  E  for  the  case  of  verti¬ 
cal  polarization  or  the  vertical  component  of  H  for  horizontal  polariza¬ 
tion.  The  time  dependence  is  exp(iu't)  and  the  normalized  impedance, 
A,  near  grazing  is 


vertical  polarization 


horizontal  polarization 

with 


o 


vhere  c:  is  the  dielectric  constant,  o  is  the  conduotivitx'  and  a  the 
r  ' 

angular  frequency. 

The  source  distribution  is  t(x, y).  Let 
cp  =  e  Ux.  y) 
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and  i)  becomes 


Assximing  that  the  fast  variation  with  x  occurs  in  exp(-ikx) 


or  tihat  d®  dst?  is  small  compared  with  remaining  terms  we  find 


^  -2TTT(x,y)  e^ 


(A-1) 


An  elementary  fxinctioa  for  (A-l)  is  (Ott  and  Berry,  1970) 


-y  —  G(x,y;§,ri)  =  - - - - 

J  %-x 


...  -ik^n 

+  — "  ■  1  exp{-ik(t-y)‘® /2(§-x)}e^  ^dt,  X  < 

./?-x  : 


,-ik(n.y)=./2(?.x) 


W(x,^).  X  < 


*  2ik 

./  G{x,  y;§,  ’I)  =  0,  x  >  %  . 

»  TT 
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The  function  satisfies 


+  2ik  =  -2tt6(x-5,  y-T))  (A-3) 

The  constant  on  the  left-hand- side  of  (A- 2)  comes  from  integrating  both 
sides  of  (A-3)  over  the  region  R  =  {  x,  y:  -  •<  x  <  y(x)  <y  <«}  . 

Multiplying  (A-1)  by  C,  (A-3)  by  ^  ,  and  subtracting  and  inte¬ 
grating  over  the  region  R  yields 


♦  17)  ♦ 


SG 

a* 


)  dx  dy 


=  -2*t 


xGdxdy+nifP) 


r 


(A-4) 


where  P  is  the  observation  point  (^.ri).  *nd  j;  is  a  region  around  the 
source.  The  divergence  theorem  yields  on  the  surface  y(x) 


^  l<ixdy=  r  (C^  - 
ay®  ay®  '  ]  "'y  =^y 

R  c 


)  In 


dc 


(A-5) 


where  e  is  the  outward  directed  normal  (into  the  surface)  and  C  is  a 
—  n 

contour  enclosing  R  and 


e 

— n 


-y*  e  +  e 

-I 

V  1  +  <y')® 
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and  along  y  s  y(x) 


dx 


Also 


2ik 


+*^dxdy=2ik  ^  ^  (C*)dxdy 
R  R 


=  2ik  j  G  'll  dc 

C 


(A -6) 


Frooi  (iiX  and  neglecting  d^/dx  compared  vddi  other  terms 

|A=  ik/ii  -  iky'  If  (A-7) 

and  substituting  (A >5),  (A >6),  and  (A«7)  into  (A>4)»  and  assuming  v  =  0, 
for  x  ^  0,  which  means  neglecting  backscatter  from  the  regiai  x  -  0  , 
and  all  sources  are  in  the  region  x  >  0  , 

^  ac  •’ 

-]  LikMG-iky''>G |^]dx-2ik  j  G  v  y' dx 
o  o 


+  2t!  jj  "  e  Gdxdy  =  *:v(P) 
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j  [-ik^i^G-iky' iiG  +  1 15  1  dx+  2n  t  G  dxdy  =  tt  MP) 


(A-8) 


Substituting  (Ott  and  Berry,  1970) 


—  —  =  ikii  /  —  G 

TT  ~  ^  ^  V  TT 


in  (A-8)  gives 


ikj  2^:1,  ax 

J  jYT::  L5.xj' 


+  2tt  rCe  d:dy  =  TT^  (P) 

i: 


•>  *“x) 

Reintroducing  "  and  defining  G  =  G  e  yields 


r;  4  - 


exp-ik((g-x)^^[(n-y)V2(S-x)  j)  n-y  ^ 


- 


+  ,  TGdxd>  -I  ^<P)c‘^ 


(A-10) 


-28- 


We  ?.s8\ixne  that  the  antenna  has  a  phase  center  v^ere  the  source 
distribution,  't(x,y),  is  located.  Then  we  write 


'^{P)  =  g(P)’^exp  + 


£»  (x,  y) 


(A-11) 


where  (x  +  y^  /  2x)  is  the  first  two  terms  in  the  binomial  expansion  for 
the  distance  between  the  source  point  O  and  the  obser^/ation  point  at  P. 
The  function  g(P)  represents  the  antenna  gain  or  pattern  facto*-.  We 
also  introduce  an  attenuation  function  f(P)  defined  as 


(I)(P)  =  Zi{h,  xp^ 


(A-12) 


When  ftese  two  equations  are  substituted  into  (A >10),  we  find  (inter - 
changing  (?,T\)  viith  (x,y)) 


^(3^  =  g(x,y)  W(x,o) 


V  2^  J  * 


\y'(?)W(x,?) 


r 

I 

u 


I 

..  c 


d 


s 


where 


*(x,?) 

y 

r 


Jbciil! 


2{x-^) 

y(x) 


2x 


(A-13) 


which  differs  slightly  from  the  result  >n  Ott  and  Berry  (1970);  sec  for 
example  Ott  (1971). 


_____  1 
yf  (x,^)  =  1  J  TT  p  s’^  erfc  (i u^) 


-  (x-§) 

2 


(A -14) 
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APPENDIX  B:  Numerical  Analysis 

The  integral  equation  (1)  or  equivalently  (A. -13)  is  of  the  form  of 
a  linear  Volterra  integral  equation  of  the  second  kindi  i.  e. , 


.X 

f{x)  =  g(x)  -  c  ]  f(s)  K(x,  s)  ds  (B-1) 

o 

■where  f(x)  is  the  imkno'wn  attenuation  function  whose  value  is  to  be 
determined  in  the  interval  0  ^  s  ^  x  .  The  function  g(x)  and  K(x,  s) 
are  kno'wn,  and  c  is  a  constant.  If  g{x)  is  bounded  and  continuous 
and  if 


J  1K(x,  s)l  ds  ^  L  <  “  (B»2) 

o 

then  the  solution  will  be  unique  and  continuous  (Wagner,  1953).  This  integral 
equation  can  be  solved  by  a  stepwise  calculation  that  divides  the  interval 
X  into  subintervals  of  arbitrary  ■width. 

That  is,  consider  the  subdivision 


.Xx 


w{x^)  -  {iA)2| 


f(s)K(x  ,  s)  ds 
n 


+  J  f{s)K  (x^,  s)dsj 


(B-3) 


The  unknown  function,  f(s),  is  fitted  with  a  polynomial  of  the  form 


f(s)  =  a^  +  ax  s  +  ag  s^ 


(B-4) 
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Increasing  the  degree  of  the  polynomial  to  3  or  higher  would 
result  in  even  higher  accuracy;  however,  the  algebra  becomes  more 
complicated  and  sufficient  accuracy  can  be  obtained  with  the  polynomial 
of  degree  2.  In  some  examples,  the  solution  may  become  unstable  for 
the  higher  degree  polynomial  and  oscillate  between  the  fitted  points. 

The  solution  of  the  integral  eqviation  requires  special  starting 
procedures.  We  suggest  that  the  interpolating  polynomial  be  of  the  form 

-  3/2 

f(s)  =  ao  +  tti  s^  +  Oa  s  +  tta  s  ,  0  ^  s  ^  Xj*  (B-5) 

and  to  use  (B-4)  for  j<c  ^  s  x^.  The  choice  in  (B-5)  is  a  logical  one  if 
we  assume  the  terrain  is  flat  in  the  immediate  vicinity  of  the  transmit¬ 
ting  antenna.  If  the  terrain  is  flat  the  exact  answer  for  the  attenuation 
function  is  then  in  fact  a  half-order  power  series  in  the  numerical  dis¬ 
tance.  Requiring  the  polynomial  in  (B-5)  to  pass  through  the  first  four 


consecutive  points  yields 

=  1.0 

(B-6a) 

O-l 

=  Ri  f(:q  )  +  Re  f(3ifc  )  +  Rg  f(x5 )  +  R, 

(B-6b) 

=  R5f(xi)  +  Re  f (3%  )  +  Rr  f(3<b  )  +  Re 

(B-6c) 

<^3 

=  Rg  f(Xi  )  +  Riof(3%  )  +  Rl  1  f(3%  )  +  RlS 

(B-6d) 

The  constants  in  (B-4]  are  found  by  requiring  the  polynomial  to  pass 
through  the  points  g  *  ^i  1  ^i  ’  ^  simple  exercise  to  show 


that 


3d  ~  Ri3  f(x.)  +  Ri4  ^is 

%  -  Ri6  f(x.)  +  Ri7  ^18 

Sr  =  Ri9f(x.)  +  Reof(x._j)  +  Rei 
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tu.j) 

(B-7a) 

(B-7b) 

(B-7c) 

where  the  R' s  in  (B-6)  and  (B-7)  are  defined  as 


if-  ii  ii  ii-i 

D  =  (Xj  3%  3%  )  ^  I^Xi  -3%“*)  +  3%  (Xi^  -3^^)  +  3%  (3%^  -Xi^)  J 

i  1 

=  XfeJfe  (313^-3%^) /D 

i  1 

Rg  =  3ei3fe  (Xi^-3{3^)/D 


Ra  =  Xi3qg  (3%^-Xi^)/D 

3/2  3/2 


„  _r  ,  3/2  3/2.,  ,  3/2  3/2.,  ,  3/2  3/2  'I,  ^ 

Rl  =LXi(3!9  -3<b  )+3{b(Xi  -3<g  )  +  3%  (sJe  -Xi  )J/D 

1 

Rg  =  (3%  3%  )  ^  (3%  -3^  )  /  D 

i 

Rg  =  (xi3%)^  (xg  -Xi)/D 

i 

R?  =  (Xl3%)^  (Xi  -3%  )/D 

„  _r  I,  3/2  3/2.  .  i,  3/2  3/2,.  i  3/2  3/2,-],_ 

Re  =|_3{;i  (3%  -3%  )+3%(Xs  -3Ci  )  +  3%‘‘(Xi  -Xfe  )J/D 

Rg  =  (3%  3%  )  ^  (3Cj^  -3%^)  /  D 
i  i  i 

R,0  =  (Xi3e,)^  (xi^-3^^)/D 

Rn  =  (xx3i;j)2  (3%^-xi2)/D 

Ria  =  -3%  )  +  3%^(3£l  “3%  )  +  3%*(3iB  -X,  )J  /  D 

D.  =  (x.  -X.  jfsj^-x.  (x.  +x.  ) 

1  1-2  1-1  L  1  1  1-2  1-1 


X.  X.  ,  , 

1-2  1-1  J 


1 


=  Vl  *i-2  ^-2 


R 


14  =  X.  X,  _  (x.  -  X.  J/D. 

1  1-2  1  1-2  1 


Ri5  =x.x._j  (x._j  -x.)/D. 
R,e  =  (=^.1 
Ri7  =^.2  )/°i 

Ris  =  (X? 


<•19  =  (*i.2 


(B-8a) 

(B-8b} 

(B-8c) 

(B  -8d) 

(B-8e) 

(B-8f) 

(B-8g) 

(B-8h) 

(B.8i) 

{B-8j) 

(B  -8k) 

(B-8f9) 

(B-8m) 

(B-8n) 

(B-80) 

(B-8p) 

(B-8q) 

(B-8r) 

(B-8s) 

(B-8t) 

{B-8u) 
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(B-8v) 

(B-8'v) 


Rgo  =  (x.  -  x.  _2)/Di 
Rei  =  (x._j  -  x.)/D. 

Using  our  polynomial  interpolation  formulas  for  f(8)  we  find 
that  the  integrals  in  (B-3)  all  have  the  following  generic  form 

(x.  ,  Xj,  =  J  s'^^^K(x.  ,  s)  ds  (B-9y 

with 

0?s  j-1 

1  ^  j  ^  i 

2  ^  i  ^  n  (B-10) 

1=0,  1,  2,  3,  4. 

These  integrals  are  evaluated  numerically  using  a  five  point  Gaussian 
integration  formula  with  special  attention  given  to  those  integrals  having 
singularities  at  either  of  the  endpoints  of  integration. 

Substituting  (B-4)  through  (B-10)  into  (B-3)  yields  the  following 
general  expression  for  f(x)  at  the  i^^  point 

f(i)  |l  +  (i/X)^[^Ri3  (i)p^(i,i,i-l)+Ri6  (i)pi3  (i,  i, i-l)+Ri9  (i)ft  (i,i,i-l) 

i  4  3  3 

=  W(i)  -{i  A  ) "  I  ^  ft,  (i,  j,  j -1)+R4  ^  Pi  (i,  j,  j  -D+Re  ^  Pi,  (i,  j,  j-i) 
j=l  j=l  j=l 

3  3  J 

+Ri2  "[r  (i,  j,  j-l)+f(l)[Ri  ^Pi  (i,  j,  j-1)  +  Ffe  (i.  j,  j-1) 

j=l  j=l  U 
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j 

i 

3  3  3 

+1^  ^Rb  (i.  j,  j-l)]+f(2)[Rfi  (i.jJ-D+Re  ^Fte  U.  j.  j-1) 

j=l  j=l  j=l 


3 

+Rio^ft  (i,  j,  j-l)  Rig  (4)Fb  (i,  4,  3)+Ria  (4)ft..  (i,  4,  3)+Rei  (4)rk  (i,  4,  3)] 

j=l 


+ 


f(3)[R3y^Pi(i,j,j-l)+R,y^Pte(i, 


+Rii 


j=l 


j=l 


j=l 


+Ri4  (4)R5  (i,  4,  3)+Ri7  (4)b2  (i,  4,  3)+Reo  (4)^  (i,  4,  3)+Ri5  (5)r5  (i,  5.  4) 


i-2 

+  Ri8  (5)R3  (i»  5,  4)+Rei  (5)i^  (i,  5,  4)  j  +  ^  Fm)  |^Ri3  (m)pb  (i, 

m=4 


+Ri4  (m+Dft,  (i,ni+l,m)  fRig  (m+2)Rj  (i, .  i+2,  iti+l) 


+Rie  (iTfi)P3  (iim,m-l)+Ri7  (m+l)£:fe  (i,m+l,m)  +Rib  (m+2)Ra  (i,m+2,m+l) 
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+  Rig  (m)pi  (i,  m, m-l)+Rog  (m+l)pj  (i,  m+1 ,  m)+RB  j  (ni+2)i^  (i,  m+2,  m+1)] 


+f(i-l)[Ri3(l-l)R5(i,i-l  1-2)  (i)j^  (i,i,i-l)+R^g  (i.l)p,  (i,i-i,i.2) 


+Ri7  (i)E3  U"^)Pfc  ^"2) +RBo(yBi(i»  (B-11) 


Reference 

(B-1)  Wagner,  Carl  (1953),  "On  the  numerical  solution  of  Volterra 
integral  equations,  "  J.  Math,  and  Phys.  32,  pp.  289-401. 
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APPENDIX  C:  The  Computer  Program  and  Flow  Chart 

Program  Wagner  implements  the  procedure  given  in  Appendix  B 
for  solving  the  integral  equation  derived  in  Appendix  A.  Flexibility  is 
obtained  by  using  appropriate  versions  of  three  subroutines: 

(1)  TERRANE,  which  returns  the  height,  slope,  and  ground  con¬ 
stants  as  a  fxinction  of  distance,  x.  By  writing  appropriate 

r 

statements  in  this  subroutine  the  user  can  define  any  propagation  path  he 
needs*  The  general  form  of  the  subroutine  TERRANE  is  sho'vn  on  page  49, 
and  two  particular  implementations  used  for  examples  in  this  report  are 
listed  on  pages  50  and  51  . 

(2)  DISTX,  vdiich  returns  the  set  of  distances  x(I)  at  which  the 
function  F(x)  will  be  calculated.  The  general  form  of  DISTX  is  shown 
on  page  45,  and  two  particular  implementations  are  shown  on  pages  46 
and  47  . 

(3)  KERNL,  vdiich  computes  the  kernel  of  the  integral  equation. 
Program  Wagner  can  be  used  to  solve  other  integral  equations  of  the 
form  (B-1)  if  the  kernel  includes  the  factor  [  8(x>s)3  by  modifying 
subroutine  KERNL.  For  example,  WAGNER  can  solve  Hufford' s 
integral  equation. 

Comment  cards  in  the  listings  that  follow  explain  the  program’s 
operation.  The  input  card  sequence  for  Prc^ram  WAGNER  is 
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Card 

Cols. 

Description 

1 

1-10 

The  number  of  Gaussian  quadrature  abcissas 
and  weights  (5  recommended) 

2  through  4 

3-33  U  36-66 

Values  for  the  Gaussian  weights  and  abcissas. 

5  through  N4-4 

1-10 

The  N  points  at  uhich  die  attenuation  func¬ 
tion  is  to  be  calculated.  These  distances 
are  read  in  IdlometerSi  by  OISTX. 

N  +  5 

A  blank  card  udiich  signals  the  end  of  the 
distance  deck  when  the  form  of  DISTX  is 
that  given  on  page  47.  When  DISTX  takes 
the  form  given  on  page  46,  no  blank  card 
ic  required. 

N  +  6 

1-10 

Source  hei^t  in  kilometers. 

11-20 

Frequency  in  Megaherts. 

21-30 

PolarisaUon,  1.  =  vertical,  2.  =  horiaontal. 

Following  is  s  flow  cliart  together  with  a  statement  listiag 
(Fortran  3800)  of  the  computer  program,  and  a  sample  output. 
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Flow  chart  for  computer  program 
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PROGRAM  WAGNER 

C 

C  A  PROGRAM  TO  COMPUTE  HF  GROUND  WAVE  ATTENUATION 
C  IRREGULAR.  INHOMOGENEOUS  TERRAIN.  REFERENCE: 

C  TELECOMMUNICATIONS  RESEARCH  REPORT.  No.  7  .  1970. 

C 

DIMENSION  IPOL(2) 

COMMON  /O/  F(2000)*R13(2000>»Rl4(2000)»NlS(2000:tR16(2000)t 
I  R17(2000)iR18(2000)*R19(2000) tR20(2300) •ft21(2000> 

COMMON  /!/  HA 
COMMON  /2/  0*H»HP 
COMMON  /3/  OELTAR.WAVE 
COMMON  tki  FREOtPOL 
COMMON  /S/  NG*AB(48)»GHU8I 
COMMON  /6/  N»X(2001)«t 
TYPE  DOUBLE  DAB *008 

COMPLEX  FEWH >F t ALAMZ (SUM tOELTARtETAR 
COMPLEX  KERNL*Po*P1*P2»P3»P4.CTMP 
1P0LI1)*8H  VERTIC  S  IPOL( 2 )«8HH0RIZ0NT 
C 

C  READ  GAUSSIAN  QUADRATURE  ABCISSaS  AHO  WEIGHTS 

C 

READ  1000*  NG 
1000  FORMAT!  no) 

NR>(NG4>ll/2 
DO  1  L«l.Nfi 
READ  1010*  DABtOGH 
1010  FORMAT! 2033*25) 

J*N6~C*1 

AB!L)-OAB 

AB!U)«>AB!L) 

GH!L)-DGM 
1  GH!J)>GM!L) 

C 

<-  CALL  subroutine  TO  SET  UP  DISTAiCE  ARRAY  X  IN  METERS 

C  START  WITH  x;2),  X!l)»0*  HAS  ALREADY  BEEN  SET. 

C  THE  DISTANCES  00  NOT  HAVE  TO  BE  E&./ALLY  SPACED. 

C  SUBROUTINE  LISTX  SHOULD  MAKE  SURE  N  <  2000 

C 

x!n*o. 

F(1)«!1.*0.) 

CALL  OISTX 
C 

C  MAKE  SURE  THEhE  ARE  AT  LEAST  A  DISTANCES 

IF  iN.GE.4)  GO  TO  4 
PRINT  1040 

1040  format  (•ONUNBER  OF  DISTANCES  0  J 
CALL  EXIT 
C 

4  S0RTX2- SORT iX! 2)1 
SORTX3*SORT!X!3n 
SQRTX4«S0RT!X(4n 

Dl*SQRTlXI2)*X!3)*X!*n«!X(2  J«!SORTX*-SURTX3»*X(3)«ISORT't2-SORTX4» 
1  ♦X!4)*lsaRTX3-SORTX2) ) 

R1»X!3»*X!4)»IS<»TX4-S0«TX3)/D1 

R2«X!2»*X!*)*!SORTX2-SQRTX4}/01 

R3*XI2»«X!3)»*SQRTK3-SORTX2J/01 

R4*!XI2)*!SU*<tX**43-SaRTX3«»3l*’Xl  3J»ISWRTX2««3-SaRTX4**3 » 

1  ♦XI4I*ISQRTX3»*3-S0RTX2»*3) )/0l 
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R5=SQRT(X(3 )*X(4) ) » ( X < 3 ) -X ( 4 ) ) /01 
R:)=SQRT(X(2i*X(4)  )*  (X(4)-X(  2  )  ) /D1 
R7=SQRT(X(2)»X(3) )* ( X ( 2 >-X ( 3 ) » /D1 

R8=(SaRTX2*(SURTX3**3-SQRTX4*«3)+SQRTX3*lSQRTX4**3-SURTX2**3) 

1  +SQRTX4»(SQRTX2**3-S0RTX3**3) )/01 
R9=S'iRT{X(3)*XU)  )♦  (  SaRTX4-SURTX3  ) /D1 
R10=SQRT«X(2)*X(4) )*(SQRTX2-SURTX4)/01 
R11=SQRT(X(2)*X(3  )  )*(SQRTX3-bORTX2 ) /D1 

R12=<SURTX2*IX(4)-X(3) ) +SURTX3* ( X ( 2 ) -X ( 4 ) ) +SQR rX4* ! X ! 3 ) -X ( 2 ) ) )/Dl 

DO  1C  M=5»N 

M1=M-1 

M2-M-2 

D2=<X(M2)-XIM1 ) )»(X(M)**2-X<M»*(X1M1)+X(M2) ) +X ( Ml ) »X ( M2 ) ) 

R13(M)=X' Ml (M2 )*(X(M2)-XiMl ) )/D2 
R14(M)=X(M)*X(M2)*(X(Mt-X(M2) )/D2 
R15 (M)=X<M)*X(M1 ) *(X(M1 )-X(M) )/D2 
Ri6(Mi=(X(Ml )**2-XtM2)**2»/02 
R17(M)=(X{M2)**2-X(M)#*2)/D2 
R18(M)=(X(M; »*2-X(Ml)**2)/02 
R19(M)  =  (X(M2)-X(M1)  »' /D2 
R20 ( M ) = ( X ( M ) -X ( M2 n / D2 
10  R21  (M)  =  (X(M1  )-X(Mn/D2 

READ  SOURCE  HEIGHT,  FRE(3UENCY,  AND  POLARIZATION 
COL  DESCRIPTION 

1-10  SOURCE  HEIGHT,  KM 

11-20  FREQUENCY,  MHZ 

21-30  POLARIZATION,  1,  *  VERTICAL,  2.  =  HORIZONTAL 

20  PEAD  2000*  HA*FREQ,POL 
2000  FORMAT  (3F10-4) 

IF  (EOF, 60  S9,22 

22  HA=HA»1,E3 
KPOL=POL 

ALAM=2»997925E2/FRE0 

WAVE=6.283ie5307/ALAM 

ALAMZ  =  (  (0«707l067812,O.To7i0678l2.'/S(JRTF(ALAM)  ) 

CALL  ,!EADING 

PRINT  2500*  FREQ, IPOL(KPOL) ,HA 

2500  FORMAT  (*OFREaUENCY  =*,F10.2 »10X,A8 »»AL  POLAR  I ZAT I ON» , iCX , *,aNTENNA 

1  height  =*,F6.2,*  METER£»// 

2  9X,*X*,14X,*Z*,  10X,»CONDUCT! VI TY* ,3X,*DI ELECTRIC* »15X,*F(X)*,22X» 

3  *TIMU4G*/8X,*(M)#,12X,*(M)*»12X»*(MhO/M)*,6X» -CONSTANT* »8X,*MAG** 

4  13X,»ARG*,16X,*(SEC)«) 

To=KLOCK(0) 

LOOP  ON  distance 

DO  100  I=1*N 

CAuL  TERRANE  (X(  I  ) ,H,HP,COND»EPS,CONDR*EPSR) 

IF  ( I.EO.l )  GO  TO  75 
D  =  Xi I )  +  (H»*2 ) /( 2.*X(  I  )  ) 

ETaR  CMPLX  (  EPSR, -17975, *C0NDR/FREQ) 

DELTAR  =  CSQRKETaR  -  1.) 
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IFilCPOL.EQ#!)  DELTAR  «  DELTAR/ETAR 
F( I )*FEWH(H*X( I ) ) 

IF  (I,LE.6)  GO  TO  75 
C 

C  J  -  2  THROUGH  A 

C 

SUM® ( 0 • *0 • ) 

DO  40  J=2»4 
P0=P1«P2®P3= (0*  *0. ) 

•C=J-1 

XP2»0.5*(X(J)+X(K)) 

XM2»0«5*(X( J)-X(K) I 
DO  35  M=1#NG 
X0=XP2+AB(M)*XM2 
CTMP=KERNL ( XO ) *GH ( M ) 

P1=P1+CTMP*SQRT(X0) 

P2=P2+CTMP*X0 
P3*P3+CTMP*SQRT ( XQ ) **3 
IF  (K.NE.l)  GO  TO  33 
X0=0»25*X( J)*{ 1,+A6(M) )*«2 
PO=PO+SQRT ( XO ) ♦KERNL ( XO ) *GH( M ) 

GO  TO  35 
33  PO=PO+CTMP 
35  CONTINUE 
P1=P1*XM2 
P2=P2»XM2 
P3®P3*XM2 

IF  (K.NE.l)  GO  TO  38 
PO*PO*SQRT(X(J) ) 

GO  TO  40 
38  P0=P0*XM2 

40  SUM=SUM+P0+R4*P1+-R8*P2+R12»P3  +F(  2  )* (  R1*P1+R5*P2+R9»P3  ) 

1  +F ( 3 ) ♦ ( R2*P1+R6*P2+R10*P3 )+F ( 4 ) * ( R3*Pl+R7*P2+Rl 1*P3 ) 

J  =  5  THROUGH  I-l 

11=1-1 

DO  50  J=5»I1 
P0=P2=P4=(0«»0.) 

XP2=0.5*(X( J)+X( J-1 ) ) 

XM2-0.5*(X( J)-X( J-1) ) 

DC  45  M=1.NG 
X0=XP2+A8(M)*XM2 
CTMP=KERNL(X0)*GH(M) 

P0=P0+CTMP 
P2=P2+CTMP*X0 
45  P4=P4+CTMP*X0**2 
PO=PO*XM2 
P2»P2*XM2 
P4=P4*XM2 

50  SUM=SUM+F ( J-2 ) ♦ ( R15 ( J ) #P0+R18 (J)*P2+R21(J) ♦PA  ) 

1  +F( J-1 )*(R14( J)*P0+R17( J)*P2+R20( J)^P4) 

2  +F(J)  *(R13( J)*P0+R16( J)»P2+R19( J)*P4J 
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J=I 

THETA=ASINF(SQRTtX{ Il)/X( !» ) ) 
CTHETA=COSF(THETA) 

P0=P2=P4=(0«»0.) 

DO  55  M=1*NG 

TEMP=1.-0.25*CTHETA**2*(1.+AB«M) )**2 
XO=X(  n»TEMP 

CTMP=SQRTCX( I ) -XO )*KERNL ( XQ ) *GH (M ) 
PO=PO+CTMP 
P2=P2+CTMP*TEMP 
55  P4=P4+CTMP*TEMP**2 

P0=P0*CThETA*SQR7(X( I ) ) 
P2=P2*CTHETA*SQRT(X( 1 ) )**3 
P4=P4*CThETA*SQRT<X( I) )»*5 


EQUATION  (B-11) 


F(  I  )  =  (F<  I  )-ALAMZ»(SUM+F(  I-2)*(R15(  I ) ♦P0-rRl8 (  n*P2+R21  ( I ) ♦P4) 

1  +F( I1)*(R14( I)*P0+R17( n*P2+R20< I )*P4) ) ) / ( 1 #+ALAMZ* ( R13 ( I) *P0 

2  +R16( I )*P2+R19< I )*P4) ) 

75  AMP  =  CABS(F{  I )  ) 

PHA  =  CANGCFd)) 

TIME=(KLOCK(0)“TO)*0.001 
PRINT  8000*  X( 1) ♦h»COND»EPS»AMP»PhA»TIME 
8000  format  (*0*»F12t2*Fl8.9»Fl4.6»F13.4>E18,8»El6.8»F15.3) 
loo  CONTINUE 
C 

GO  TO  20 
999  CALL  EXIT 
END 
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SUBROUTINE  DISTX 

C  READ  DISTANCES  IN  KM  AND  CONVERTS  THEM  TO  METERS 

C  (A  DISTANCE  OF  ZERO  SIGNALS  END  OF  DISTANCE  DECK) 

COMMON  /6/  N»X(2001)*I 


IN  THIS  SUBROUTINE  THE  USER  MUST  FILL 
THE  X(I)  ARRAY  WITH  N  VARIABLES. 


RETURN 

END 
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SUBROUTINE  OISTX 

C  COMPUTES  EQUALLY  SPACED  DISTANCES 

COMMON  /6/  N*X(2001J»I 
C  INPUT 

C  DMIN  —  FIRST  DISTANCE  IN  KM 

C  DMAX  — •  MAXIMUM  DISTANCE  IN  KM 

C  DINC  INCREMENT  ON  DISTANCE  IN  KM 

C 

READ  1000*  DMIN»DMAX»DINC 
1000  FORMAT  (3F10.2) 

IF  (DMIN.EQ*0*)  0MINBDMIN4-DINC 

N=(DMAX-DMlN)/DINC+2 

00  10  I»2*N 

X( I )=(OMIN+( I“2)*0INC)»1*E3 

10  continue 
RETURN 
END 


Note,  this  is  an  example  of  subroutine  DISTX 


r»  r» 


SUBROUTINE  DISTX 

READ  DISTANCES  IN  KM  AND  CONVERTS  THEM  TO  METERS 
(A  DISTANCE  OF  ZERO  SIGNALS  END  OF  DISTANCE  DECK) 
COMMON  /6/  N*X(2001)*I 
DO  2  1=2*2001 
READ  1020*  X( I ) 

1020  format  (F10»5) 

IF  (X( I ).EQ.O*)  GO  TO  3 
X(  I  )=X(  n*l*E+3 

2  CONTINUE 
PRINT  1030 

1030  FORMAT  (*0NUMBER  OF  DISTANCES  EXCEEDS  DIMENSION*) 

CALL  EXIT 

3  N=I-1 

END 


Note,  this  is  an  example  of  subroutine  DISTX, 
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FUNCTION  KERNL(XO) 

C 

C  SUBROUTINE  OF  WAGNER.  COMPUTES 

C  KERNEL  OF  INTEGRAL  EQUATION.  SEE 

C  EQ.  (A-13). 

C 

COMMON  /!/  HA 
COMMON  11/  0»H»HP 
COMMON  /3/  DELTAR.WAVE 
COMMON  /4/  FREQ»POL 
COMMON  /5/  NGiAB(48) ♦GH(48) 

COMMON  /6/  NX(iX(200l )  ♦  I 

COMPLEX  KERNL»FEWH»0ELTA*0ELTAR»ETA 

CALL  TERRANE(XO»HO»HPO»COND»EPS*CONDR»EPSR) 

ETA»CMPLX( EPS. -17975. *COND/FREQ) 

DELTA*CSQRT(ETA-1.) 

IFiPOL.EQ.l. )  DELTA=DELTA/ETA 

XMS*X( I )-X0 

HO*H-HO 

Rrt  «^WAVE*(X0  (H0**2>/ (2**X0 I )  +  XMS  +  ( (H0**2) / ( 2.»XMS) >  “  D) 
KERNL»CMPLX ( COSF «  RW ) .-SINF ( RW ) ) *SQRT ( X ( 1) / ( R1*XMS ) ) ♦ ( ( HPO+DELTA 
1  -DELTAR)*FEWH(H0»XMS»  -  (HD/XMS)) 

C 

C  THE  FACTOR  (DELTA -DELTAR)  ARISES  IN 

C  MIXED-PATH  PROBLEMS. 

C 

RETURN 

END 


SUBROUTINE  TERRANE  (X»H*HP*CONO*EPS*CONDR»EPSR) 

C  SUBROUTINE  FOR  WAGNER.  DEFINES  TERRAIN,  PROFILE  AND 

GROUND  constants. 

C 

C  INPUT  IS  DISTANCE  X  IN  METERS. 

C  OUTPUT  IS  TERRAIN  HEIGHT,  H,  SLOPE,  HP, 

C  GROUND  CONSTANTS,  CONDR,  EPSR,  COND,  EPS. 

C  IN  MIXED  PATH  CALCULATIONS,  CONOR  AND  EPSR 

C  ARE  relative  VALUES  FOR  0  AND  e  . 

C  THEY  ARE  USED  TO  COMPUTE  ' 

C  DELTAR  IN  PROGRAM  WAGNER. 

C  IN  FUNCTION  KERNL  THE  DIFFERENCE 
C  (DELTA -DELTAR)  IS 

C  COMPUTED,  this  DIFFERENCE  TAKES  INTO 

C  ACCOUNT  CHANGES 

C  IN  a  AND  Sy  WITH  DISTANCE. 

C  CONDR  AND  EPSR  ARE  USUALLY 
C  TAKEN  TO  BE  THE  VALUES  OF 

C  a  AND  e  FOR  THE  FIRST 

C  SECTION  OF  PATH. 

C 


C 

C  IN  THIS  SUBROUTINE  THE  USER  MUST  DEFINE  THE  FOLLOWING 
C  VARIABLES 

C  H  = 

C  HP  = 

C  CONDR  = 

C  EPSR  = 

C  COND  = 

C  EPS  = 

C 

C  PRINT  HEADING 

ENTRY  HEADING 
PRINT  50»A 

50  FORMAT  (*A  SMOOTH  SPHERE  WITH  RA0IUS*»E12.3) 

RETURN 

END 
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subroutine  TeRRANE  (X*H*hP*COND*ERS»CONDR*EPSR) 

SUBROUTINE  FOR  WAGNER*  DEFINES  TERRAIN. 

SMOOTH  SPHERE 
COMMON  /!/  HA 
DATA  (A-8.5E6) 

COMPUTE  HEIGHT .SLOPE .CONDUCTIVITY  AND  DIELECTRIC  CONSTANT  AT  X 
HP«-X/A 
H«.5*X*HP-  HA 
CONOR  ■  .01 
EPSR  *  10. 

CONO  »  .01 
EPS  »  10. 

RETURN 

PRINT  HEADING 
ENTRY  HEADING 
PRINT  50. A 

FORMAT  (*A  SMOOTH  SPHERE  WITH  RAD IUS*»E12 .3 ) 

RETURN 

END 


Note,  this  is  an  example  of  subroutine  TERRANE 


r»o  onr>  rior>r> 


subroutine  TERRANE  (X*H»HPtCON0»EPS*CONDR»EPSR) 

SUBROUTINE  FOR  WAGNER*  DEFINES  TERRAIN* 

TABLE  mountain  PATH  WITH  KBOL  AS  TRANSMITTER 
COMMON  /!/  HA 

compute  MEIGMT»SL0PE»C0N0UCTIVITY  and  DIELECT'^IC  constant  AT  X 
H  «  50**TANHF((X-5000*>/100*»+50.-HA 
HP- 0 . 5* ( 1 * - 1  TANhF ( ( X-5000 . ) / 1 00  * ) » **2 ) 

CONOR  -  *01 
EPSR  -10* 

A  FOUR  SECTION  PATH 

XI  «  28574*0 
X2  ■  35000* 

X3  s  45000*0 

1F(X*GT*X1*AN0*X*LE*X2)  60  TO  20 
IF(X*6T.X2*AND*X*LE*X3)  GO  TO  30 
IF(X*6T*X3)  GO  TO  40 
CONO  »  .01 
EPS  >  10* 

GO  TO  10 
20  CONO  «  2.0 
EPS  ■  .81 
GO  TO  10 
30  CONO  «  *01 
EPS  ■  10* 

GO  TO  10 
40  CONO  -2.0 
EPS  «  81.0 

10  CONTINUE 
RETURN 

PRINT  HEADING 
ENTRY  HEADING 

PRINT  50  „  . 

50  FORMaT(»TACLE  MOUTAIN  Path  with  XBOL  aS  transmitter*! 

RETURN 

END 


Note,  this  is  an  example  of  subroutine  TERRANE. 
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COMPLEX  FUNCTION  FEWH^HOtXD) 

C 

C  THE  ATTENUATION  FUNCTION, 

C  EQ  (A-U),  OF  TELECOMMUNICATIONS  RESEARCH 

C  REPORT  No.  7  .  1970.  INPUT  IS  THE 

C  HEIGHT  HD  AND  THE  DISTANCE 

C  XD. 

C 

COMMON  /3/  OELTARtteAVE 

COMPLEX  FEWH  *  TEMP  *0  tltlZ tZZ  tHWERF • WERFZ • WERF  tZVERF  tDELT AR 
TEMP»( 0. 7071067812 f-0* 7071067812 )*SORT<.5*WAVE) 

X02«SQRT(XD) 

Q»-TEMP*H0/X02 
Z«TEMP*DELTAR*X02  ♦  0 
Z2— Z 

Z!«AIMAG(ZZ) 

IF  (ZI.LT«0.«UR.(A6S(REAL(ZZ)  )«LT«6..AND*ZULT.6.n  60  TO  lo 
Z2»ZZ**2 

HWERF»(Z2-2.»/UZ»(Z2-3*5)  ) 

GO  TO  12 

10  yERFZ»WERFlZZ) 

HWERF«ZZ-0*5*WERF2/<ZZ*WERF2+(0.»-0*56418958» ) 

12  ZWERF-Z+HWERF 

FE»(H»  (Q*ZWERF~0.5)/(Z»2WERF-0.5  > 

return 

END 
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COMPLEX  function  WERFIZZZ) 

C 

C  THE  FUNCTION  w(7,), 

C  ABRAMOWITZ  AND  STEGXB^,  1964) 

C  WRITTEN  BY  DR.  GEORGE  HUFFORD,  AND  MODIFIED  BY 
c  dr.  R.  H.  OTT 

C 

COMPLEX  Z*ZZZ»ZV*V*Z2»CtW*S 
DIMENSION  C(12)»W(5«4) 

equivalence  <S»C<12n 

DATA  CCd)  «  (•0»'‘*664189S83S)  ) 
data  (( (W( I »J) • I«1*3) • 

X  (3.6787944ll714423E-0lt6.071577068413937E“01 ) . 

X  ii.83iS638888734l8E-02*3*400262170660662E-01 ) . 

X  {1.2340v80‘»0866788E-04.2.011673170376004E-01)  . 

X  (1.125351747192646E-07.1.459535899001528E-01 ) * 

X  (4.27S83S76L3S8070£-01*0*OOOOOOOOOOOOOOOE+00) » 

X  (3.047442032569i26£-01.2.082189382028316E-01) . 

X  { 1.4023958 13662  n9E~0l  *2. 222 134401798991E--01 )  . 

X  (6.531777728904697E-02*1.739183154163490E-01  )  * 

X  (3.628145648998864E*02.1.3583e95i0006551E-01  )  . 

X  (  2 .553‘'  6763105058E-01  •O.COOOOOOOOOOOOOOE+00  )  t 
X  (2.1849261527489076-01*9. 299780939260186E-02 ) » 

X  (1.47952  55120158E-01*1*31179717o842178E-01) » 

X  (9. 271076*  •*26443326-02. 1. 283169622282615E-01 )  ♦ 

X  (5.96869296lr44590£-02.1.132l00561244882E-01) . 

X  ( 1.7900li5il8i3930E-OltO.OOOOOOO-30000000E'^00)  ♦ 

X  (  1.6‘*26li36i92986l£-0l.5.0l97135l3524966E-02J  . 

X  (1.307574696698522E-01.8.111265047745472E-02 ) * 

X  (9.640250558304439E-02.9.123632600421258E-02 ) » 

X  (6.979o96164964750£-02.8.934000024036461E-02 ) ) 
XX»REAL{ZZZ) 

YY=AIMAG(ZZZ) 

X«ABSF(XX) 

Y=ABSF(YYI 

Z»CMPLX(X.Y} 

LZ2=0 

1F(X.GE.4.5.0R.Y.OE.3.5)  GO  TO  loO 


CONVERGING  SERIES 
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I=X+.5 

J=Y+.5 

V*CMPLX(FL0ATF( I) ♦FLOATCJJ  J 
ZV=Z-V 

C(2)«W(  I-t-l.J+1) 

AI*0» 

DO  10  1*3»12 
AI=AI>.& 

C(  n  =  (V*C(  I-l)+C{  1-2)  )/AI 

10  CONTINUE 
J=12 

DO  11  1=2*11 
J=J-1 

11  S=S*ZV+C(J) 

20  IF(YYcGE.0.>  GO  TO  30 
IF(.NOT.LZ2)  22=Z*Z 
S=2.*CEXP(-Z2)-S 
IFiXX.GT.O.l  S=C0NJ6(S) 

GO  TO  200 

30  IF<XX.LT,0.)  S=CONJG(S) 

200  WERF=S 
RETURN 
loo  L22=l 
Z2*Z»Z 


C 

C  ASYMPTOTIC  SERIES 
C 


S  =  Z»((0..0.4613135279)/(Z2  -  0.1901635092)  +  ( 
X(Z2  -  1.7844927485)  +  < 0. *0.0028838938748 )/( 22  - 
GO  TO  20 
END 


0. *0.09999216168)/ 
5.52534374379)  ) 
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Input  data  for  the  case  of  a  smooth  cylindrical  earth; 

Card  #1:  5  (Column  5) 

Card  #2:  0.9061798459  0.2369268851 

Card  #3:  0.5384693101  0.4786286704 

Card  #4:  0.0000000000  0.  56888888888 

Card  #5 

through  62:  1.0,  2.0,  . .  53. 

Card  #63:  0.0  (column  8)  1.0  ^column  18)  1.0  (column  28) 

Following  is  the  output  from  this  example. 


55 


A  SHOOTh  SPHEKt  WITH 


S.SOQ'taOo 


frequenct  = 

1.  0  J 

X 

2 

(M) 

(M) 

0.00 

0. ooooooaoo 

1000.00 

-0.058823529 

2000. 00 

-0.2352 94110 

3000.00 

-0.529411765 

4000.00 

-0.9+1176471 

5000.00 

-1.4705o8r35 

6000.00 

-2.  117647059 

7000.00 

-2. 352352941 

8000.00 

•3. 734705802 

9000.00 

-4.  734705082 

10000.00 

-5.852352941 

1  lOOC.OO 

-<•.  117647059 

12000.00 

•8.470580235 

13000.00 

-9.9+11 76471 

14000.00 

-11.529411765 

15000.00 

-13. 2T52  94118 

16000.00 

-15.  0  >8823529 

17000.00 

-1  7.  0  000  000  00 

18000.00 

-19.0  18023530 

19000.00 

-21.2 15294118 

20000.  00 

-23. 579411765 

21000.00 

-25.9+1176471 

22000.00 

-2b.  +70586239 

2  3  u  j  0  ,  0  0 

-  71.11  7647058 

24000.00 

-33. 852352941 

25C0i:.  00 

-36.7  <470  188, 7 

VERTICAL  POLARIZATION 


conductivity 

IMHC/M) 

dielectric 

CONSTANT 

8  .019008 

1  0.  0080 

0.010000 

10.8000 

8,01000  0 

10.0000 

0 .010000 

10.0000 

0.010000 

in.oooo 

O.OIDOOO 

10-0000 

O.OICOOO 

10.0000 

0  .010000 

10.0000 

0 .010000 

10.0000 

0.010000 

10.0000 

0 .010009 

10.0000 

0.010800 

10.0000 

9.018000 

10. 0000 

0.010000 

1  0.  0000 

0.010300 

10.  0000 

8.010000 

1 0.0000 

0.010000 

10.0000 

Q.CIOOOO 

1  0.0000 

0.010000 

1  0.  0000 

0.010000 

10. 0000 

a  .01  0008 

10. 0000 

O.ClOOOO 

10. 0000 

0  .DtOO'll 

1 0. oooo 

8.010000 

10.0000 

0.010803  .5^. 

10. 0000 

ANTENNA  HEIUHT  =  0.00  METERS 

FIX) 

MAO  ARG 

l.OOODOOOOt-llOO  O.OOOOOOOOvOOO 
9.627S5706-001  -4.24608707-001 
9.34092461-001  -5.97768976-001 
9.07383363-001  -7.28975441-001 
8.81991252-001  -8.38214925-UOl 
8.57643738-001  -9.33250508-001 
8.34825345-001  -1.02254934+000 
8.12513417-001  -1.10098530+000 
7.90993050-001  -1.17332944+000 
7.70208750-001  -1.24066203+000 
7.50115518-001  -1.30377135+000 
7.30675486-001  -1.36325394+000 
7.11855902-001  -1.41957439+000 
6.93627838-001  -1  .47310306  +  000 
6.75965320-001  -1.52414101+000 
6.58844727-001  -1.57293697+000 
6.42?4'.356-001  -1.61969937  +  000 
6.26144100-001  -1.66460498+000 
6.10525209-001  -1.70780530+000 
5.95370101-001  -1 . 74943145+0 00 
5.80662214-001  -1.78959780+000 
5.66385892-001  -1.82840489+000 
5.52526281-001  -1.86594171+000 
5.39069260-001  -1.90228746+000 
5.26001367-001  -1.93751312+000 
5.13309747-001  -1.97168246+000 


timing 

(SEO 
0.000 
0.007 
0.033 
0.045 
0.056 
0.069 
0.348 
0.665 
1.022 
1.413 
1.838 
2.310 
2.819 
3.358 
3.942 
4.558 
5.216 
5.913 
b  •g4»9 

6*231 
9.094 
9,982 
1 0  •  90  7 
-i.au9 
12,674 


26000.  .0 

-39.7j47U58#2 

8 .Ui U8UU 

lUaUUUU 

2  7000.00 

-42.312352942 

0  .01  0000 

10.0000' 

28000.00 

-46.117647058 

0  .01  0000 

10.0000 

29000.00 

-49.470588235 

0.01  0000 

10.0000 

30000.00 

-52.9*1176470 

0  .01  0300 

10.8000 

31000.00 

-56.529411764 

0.010000 

10.0000 

32000.00 

-60.235294118 

0  .010000 

10.0000 

33000.00 

-64.055823530 

0 .01000  0 

10.  0000 

34000.00 

-68.030000000 

0  .01000  0 

10.  0000 

35000.00 

-72.058823530 

0.010000 

10.0000 

36000.00 

-76.235294119 

0  .01000  0 

10.0000 

37000.00 

-80.529411763 

0.010000 

10.0000 

38000.00 

-84.9*1176470 

0.010000 

10.  0000 

39000.00 

-89.470588237 

0  .010000 

10.0000 

40000.00 

-94.  117647059 

0  .010000 

10.0000 

41000.00 

-90.882352941 

O.CIOOQQ 

10.0000 

42000.00 

-103.754705081 

0  .010000 

10.0000 

43000.00 

-108.764705881 

0.010000 

10.0000 

44000.00 

-113.882352941 

0.010000 

10.0000 

45000. 00 

-119.117647059 

0.013000 

10.0000 

46000.00 

•124.470588233 

0.010000 

10.0000 

47000.00 

-129.941176470 

O.CIOOOO 

10.0000 

48000.00 

-135.529411763 

0.010000 

10.0300 

49000.00 

-141.235294119 

0.01  003  0 

10.0000 

50000.00 

-147.058023530 

0.010000 

10.0000 

•  COj.OC 

i''.ooooaoooa 

o.oinj": 

10.0000 

52000.00 

-159.058823530 

O.CIOOOO 

1  0.  0030 

53000.00 

-165.235294115 

0.810000 

10.  0000 

54000.00 

-171.529411767 

O.CIOOOO 

10.0000 

55000.00 

-177.941176470 

0  .010000 

10.0000 

56000.00 

-184.4785  88233 

O.CIOOOO 

10.0000 

57000.00 

-191.  117647056 

0,01  0303 

10.0000 

58000.00 

-197.802352941 

0.010000  -57- 

10.0000 

9  jui  -i.  j»u  uu  la.itw 

‘•4a90066lt9*D0l  -2.03707?62»000  15.002 

4.27372092-001  -2. 0 6040402*0 00  16.113 

4.66067577-001  -2.05607556*000  17.270 

4.55002679-001  -2.12053262*000  16.962 

4.44407365-001  -2.15741220*000  19.696 

4.34031904-001  -2.10554055*000  20.971 

4.23947239-001  -2.21297341*000  22.274 

4.14144170-001  -2.23971633*000  23.624 

4.04614161-001  -2.26560409*000  25.009 

3.95340947-001  -2.29126493*000  26.429 

3.66340400-001  -2.31612073*000  27.091 

3.77579473-001  -2.34030399*000  29.399 

3.69060501-001  -2.3  640  9375*0  00  30.941 

3.60775570-001  -2.30726194*000  32.510 

3.52717656-001  -2.40990741*000  34.150 

3.44079922*001  -2.43204799*0  00  35.003 

3.37255794-001  -2.45370059*000  37.490 

3.29030903-001  -2.47400123*000  39.216 

3.22623091-001  -2.49560513*000  40.990 

3.15602395-001  -2.51500601*000  42.606 

3.08771045-001  -2.53574010*000  44.655 

3.02123453-001  -2.55517017*000  46.549 

2.95654207-001  -2.57421367*000  40.476 

2.09350067-001  -2.59205067*000  50.431 

2.03229953-001  -7.61112474*000  52.444 

2.77264947-001  -2.62902301*000  54.401 

2.71450201-001  -2.64656414*000  56.555 

2.65005333-001  -2.66375042*000  50.664 

2.60301624-001  -2.60061575*000  60.010 

2.54942015-001  -2.69714566*000  63.016 

2.49724697-001  -2.71335735*000  65.245 

2.44643190-001  -2.72925972*000  67.516 


